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We present a spin-rotation-invariant Green-function theory of long- and short-range order in 
the ferro- and antiferromagnetic Heisenberg model with arbitrary spin quantum number S on a 
stacked square lattice. The thermodynamic quantities (Curie temperature Tc, Neel temperature 
Tjv, specific heat Cv, intralayer and interlayer correlation lengths) are calculated, where the effects 
of the interlayer coupling and the 5 dependence are explored. In addition, exact diagonalizations 
on finite two-dimensional (2D) lattices with S ^ 1 are performed, and a very good agreement 
between the results of both approaches is found. For the quasi-2D and isotropic 3D magnets, 
our theory agrees well with available quantum Monte Carlo and high-temperature series-expansion 
data. Comparing the quasi-2D S = 1/2 magnets, we obtain the inequalities Tjv > Tc and, for small 
enough interlayer couplings, Tjv < Tc- The results for Cv and the intralayer correlation length 
are compared to experiments on the quasi-2D antiferromagnets Zn2VO(P04)2 with S' = 1/2 and 
La2Ni04 with S = 1, respectively. 
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I. INTRODUCTION 

Low-dimensional ferromagnetic (FM) and antiferro- 
magnetic (AF) quantum spin systems,— such as the quasi- 
two-dimensional (2D) Heisenberg ferromagnets [e.g., 
K2CUF4 with spin S = 1/2 (Ref. 0)] and antiferro- 
magnets [e.g., La2Ni04 with spin 5* = 1 (Ref. 0) being 
isostructural to the high-Tc parent compound La2Cu04] , 
are of current interest. Their study is motivated by the 
progress in the synthesis of new low-dimensional mate- 
rials. For example, very recently a defective graphene 
sheet was reported to be a room-temperature ferromag- 
netic semiconductor that may be described by an effec- 
tive quasi-2D Heisenberg model4 

Investigations of layered Heisenberg magnets by nu- 
merical methods, e.g., quantum Monte Carlo (QMC) 
simulations and high-temperature series expansions (SE), 
have been performed for a selected number of cases and 
quantities only. QMC data are available for quasi-2D 
and spatially isotropic 3D antiferromagnets with S — 1/2 
(Refs.iiH) and S' = 1 (Ref.i). SE results exist for the 
3D antiferromagnet with S = 1/2, 1, and 3/2 (Ref. §) 
and for the 3D ferromagnet with = 1/2 (Refs.[^ and[0) 
and S = 1 and 3/2 (Ref. |8|). Note that numerical stud- 
ies of ferromagnets and of 5 > 1/2 systems are rather 
scarce. 

On the other hand, analytical approaches which are 
capable to evaluate the thermodynamics of layered ferro- 
and antiferromagnets with arbitrary spin below and 
above the magnetic transition temperature Tm [M — 
C, N; Tc {Tn) denotes the Curie (Neel) temperature 
in the FM (AF) case] are desirable. In particular, the 
relation between Tm and the relevant exchange cou- 
plings can be used to determine those couplings from 
experiments. Moreover, analytical theories may have 



the advantage of being applicable in such cases, where 
the QMC method cannot be applied, e.g., in the pres- 
ence of frustration. However, the mean-field spin- 
wave theories based on the random-phase approximation 
(RPA)fiSiii that is equivalent to the Tyablikov decou- 
pling of Green functionsji^ and on auxiliary-field rep- 
resentations (Schwinger-bosoufi^iii Dyson-Maleev^i^ and 
boson-fermion representations^^) are valid only at suffi- 
ciently low temperatures and do not adequately take into 
account the temperature dependence of magnetic short- 
range order (SRO) in the paramagnetic phase. For the 
3D antiferromagnet, this deficiency has been removed by 
the quantum hierarchical reference theory of Ref. [l3- For 
quasi-2D ferro- and antiferromagnets, an essential im- 
provement in comparison to the standard mean-field ap- 
proaches may be achieved by employing the second-order 
Green-function technique^^ that we call, in the absence 
of spin anisotropics, rotation- invariant Green- function 
method (RGM). This technique provides a good descrip- 
tion of SRO and long-range order (LRO) and has been 
applied recently successfully to low-dimensional quantum 
spin ^yctcmc .-^^-^'-'-^^-^^-^'^-^^-^^-^^'^''-^^'^^-'^" 

In this paper we use the RGM and develop a theory 
of magnetic order in ferro- and antiferromagnets on a 
stacked square lattice. Thereby, we extend the previous 
work on the quasi-2D 5=1/2 antiferromagnet^^ and the 
layered S* = 1/2 ferromagnet^ to arbitrary values of the 
spin quantum number. We perform a systematic study 
of thermodynamic properties, where we contrast the FM 
with the AF cases. This allows to explore the role of 
quantum fluctuations. 

We consider the 3D spatially anisotropic Heisenberg 
model with arbitrary spin S, 



J|| ^ 
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[{'iij)xy and {i,j)z denote nearest-neighbor (NN) sites in 
the xy plane and along the z direction of a simple cubic 
lattice, respectively] with Sf = S{S + 1). For the layered 
ferromagnet (antiferromagnet) we have < (J^ > 0), 
where ^ —\\, _L. We calculate the thermodynamic prop- 
erties (magnetic transition temperatures, specific heat, 
and correlation lengths) and study the crossover from 
isotropic 2D {J± ~ 0) to 3D {J± = J||) quantum mag- 
nets. For comparison, we perform Lanczos exact diag- 
onalizations (ED) to calculate the ground state of the 
2D antiferromagnet with S* = 1, |, and 2 on a lattice 
of iV = 16 sites and full ED to get the thermodynamic 
quantities for the 2D S — 1 ferromagnet on a lattice of 
TV = 8 sites. 

The rest of the paper is organized as follows: In Sec, [ill 
the theory based on the RGM for model JT]) is developed, 
where the extension of previous RGM approache o^^i^^ 
to arbitrary spins implies novel technical aspects. In 
Sec. mil the thermodynamic properties of the 2D and 
3D ferromagnets and antiferromagnets are investigated 
as functions of temperature, spin, and interlayer cou- 
pling, also in comparison to available QMC and SE data, 
and are related to experiments. Finally, a summary of 
our work is given in Sec. IIVI 

II. ROTATION-INVARIANT 
GREEN-FUNCTION THEORY 

To evaluate the spin-correlation functions and the 
thermodynamic quantities, we calculate the dynamic 
spin susceptibility Xq~{^) ~ ^{{Sq',SZq))Lu (here, 
{{...]... ))^ denotes the two-time commutator Green 
function^^) by the RGM,^* Using the equations of motion 
up to the second step and supposing rotational symmetry 
in spin space, i.e., (Sf) = 0, we obtain uj^{{Sq; SZq))Lu = 

Mq + {{-S+;SZg))^ with Mq {[[S+,HlSZg\) and 
-S+ ^ [[S+,H],H]. For the model ^ the moment Mq 
is given by the exact expression 

Mq = -8J||Cioo(l -7q) -4J±Cooi(l -cosg,), (2) 

where Cmni = Cr = (S+S;^) = 2(5g5|j), R ^ me, + 
ney + lsz, and 7q = ^{cosQx + cosqy). The second deriva- 
tive —S^ is approximated in the spirit of the schemes 
employed in Refs. HlHIiilii, andlM That means, in 
—S^ we decouple the products of three spin operators 
along NN sequences as 

S+S+S[- = ai^{S+Si)S+ + a2^{S+Si)S+, (3) 

where the vertex parameters ai^ and Q!2^ are attached 
to NN and further-distant correlation functions, respec- 
tively, either within a layer (/i =||) or between two lay- 
ers (/i =-L). The products of three spin operators with 
two coinciding sites, appearing for 5 ^ 1, are decoupled 
asiS^^ 

S+Srs+ = {Srs+)S+ + X^{S+S-)S+, (4) 



where the vertex parameter is associated with the NN 
correlator in the layer or between NN layers. We obtain 
-S+ = u;lS+ and 

xJ-(^) = -{{S+;SZ,))^ - ^j^^, (5) 

with 

= (l-7g){A|| +16J|2ai||Cioo(l-7g)} 

+ {1 - coaqz){A± + 4j]^ai_LCooi(l - cosq^)} 
+A(l-7g)(l-cosg,), (6) 

A|| = 2j2{5 + 2A||Cioo + 2a2||(2Ciio + C2oo) (7) 
- 10q!i||Cioo} + 8J|| Jj_(q!2±Cioi — ai||Cioo), 

= Ji{^ + 2AiCooi + 2a2iCoo2 -6aiiCooi} 
+ 8J|| J_i(q!2±Cioi - ai^Cooi), (8) 

A = 8J||Ji(ai||Cioo + ai_LCooi), (9) 

where S — ^S{S + 1). From the Green function ([5]) the 
correlation functions Cr = jj Cqe^"^^ are determined 
by the spectral theorem,— 

Cq = {S+SZ,) = ^Jl + 2niu;q)], (10) 

where n{uj) — {e^^"^ — is the Bose function. The NN 
correlators are directly related to the internal energy u 
per site, u = 3J||Cioo + § Jj.C'ooi, from which the specific 
heat Cy = du/dT may be calculated. Taking the on-site 
correlator Cr=o and using the operator identity Sf — 
S+S- - + {Sf f, we get the sum rule 

Let us consider the static spin susceptibility Xq = Xq{^ = 

0) with Xq(w) = XqH^) = i-e-- Xq = Mq/2t^q- 

The lowest-order expansion of Mq and uj^ at q — yields 
Xq = HqI + ql) + hql]/[c{ql + g^) + dql], where a = 
-J||Cioo, h = -JxCooi, c = A||/4, and d = A±/2. Cal- 
culating the uniform static susceptibility x = limq_,o Xgj 
the ratio of the anisotropic functions Mq and must 
be isotropic in the limit q ^ 0, i.e., limg^^^j^o Xglq^^o — 
limq^^o Xq\q^{y)=o- That is, the condition a/c = b/d has 
to be fulfilled which reads as the isotropy condition 

4 2 

X = --^•^yC'ioo = --r — J±Cooi- (12) 
A|| A_L 

Note that such a condition was also employed in Refs.[20l. 
andii 
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The phase with magnetic LRO at T ^ Tm is described 
by the divergence of the static susceptibihty at the or- 
dering vector go, i-e-, by Xqo ~ ^'^^^ 90 = and 
Qo = Q = ('"■, I") in the FM and AF case, respectively. 
In this phase the correlation function Cr is written aai^ 



(13) 



q(#<jo) 



with Cq given by Eq. (fTO|) . The condensation part C de- 
termines the magnetization m that is defined in the spin- 
rotation-invariant form 7v? 



2Af 



The LRO conditions for the ferromagnet and antiferro- 
magnet read as = [cf. Eq. p^ ] and wq — 0, 
respectively. 

The magnetic correlation lengths above Tm niay be 
calculated by expanding Xq in the neighborhood of the 
vector go-— "^"^ For the ferromagnet (go — 0), the ex- 
pansion yields Xq = x[l + (^^ + 9y) + ^i?^]"^ with the 
squared intralayer =||) and interlayer (/i =\-) correla- 
tion lengths 



(14) 



For the antiferromagnet, the expansion around Qq = Q 
gives Xq = Xq[1 + ^likl + fc^) with k = q Q 

and 



1 



S^iiC'ioo 



-^(A||-f644ai||Cioo + 2A)-^^Ji^, (15) 



Q 



^ (A^ + 16Jiai^Cooi+2A)-^4^^- (16) 



Q 



To evaluate the thermodynamic properties, the corre- 
lation functions Cr and the vertex parameters ai^, a2^, 
and appearing in the spectrum ujq [Eqs. (H])-®] as well 
as the condensation term C in the LRO phase have to be 
determined. Besides Eqs. pU)) and (fT^ for calculating 
the correlators, we have the sum rule (fTTjl . the isotropy 
condition (H^), and the LRO conditions for determining 
the parameters; that is, we have more parameters than 
equations. To obtain a closed system of self-consistency 
equations, we reduce the number of parameters by rea- 
sonable simplifications that we have to specify for the 
FM and AF cases. 

(i) Ferromagnet: Considering the ground state {T — 0), 
we have the exact result 



Cr{0) = r/5R^o + 



(17) 



which can be reproduced by Eq. (fT5|) . Cr{0) = 
iEqi^o)[MqiO)/2iOqme^i^ + CiO), if C(0) = 



N 

and Mq{0)/2ujq{0) = jS. The equality M^{0) = 
^S^ujq{0) requires the equations ai^(O) = | and 
A;^(0) = (LRO condition, see above) or, explicitly. 



J|l (l + 1 + A|| + 3a2 



15^ 
2 



2J± {a2i_ 



J± {l + ^ + \± + a2± - I) + 4J|| (a2_L - i) = 0. In the 
special case 5 ~ 1/2, in —S^, products of spin oper- 
ators with two coinciding sites do not appear, which 
is equivalent to setting = 0. Then, the solution 
of the equations A^(0) = yields (0) — |, i.e., 
we have a2/i(0) = q;i^(0). We take this equality also 
for 5^1 and get Xfj_{0) = 2 — ^. To determine 
the parameters at finite temperatures, we first consider 
the high- temperature limit, where all a parameters ap- 
proach unity,^^ limr^oo CKi,2ft(T') — 1, and the high- 
temperature series expansior>2Z yields limx^oo Xfj.{T) = 
Aoo = 1 — 3[45(5 + 1)]"^. Because we have identical 
vertex parameters a2p, and as well as identical pa- 
rameters A|| and Aj_ at T = and for T —>■ 00, we put 
02^(7^) = aiM(T^) = ^t^iT) and A||(r) = X^{T) = A(T) 
in the whole temperature region. Then, at T ^ Tc we 
have the four parameters aj_, A and C. For their 
determination, besides the sum rule pT|) and the LRO 
conditions, A|| = and Aj_ = 0, we need an additional 
condition. Reasoning similarly as in Ref. [l^ for a param- 
eters, we consider the ratio 



A(r) - Aoo 
a||(r)-l 



^a(O) 



(18) 



as temperature independent. For T > Tc (C — 0) we 
have A^ > 0, and the number of quantities and equa- 
tions [Eqs. (HU, (dH), (UHl)] is reduced by one. 
(ii) Antiferromagnet: As revealed by previous studies of 
the 2D 5=1/2 antiferromagnet contrary to the FM 
case, the introduction of the vertex parameter a2 7^ ai 
appreciably improves the results as compared with the 
simplification a2 — cti. We expect the same behavior 
also for the layered antiferromagnet. This can be under- 
stood as follows. In the LRO phase and paraphase with 
AF SRO, the parameter ai^ is associated with NN cor- 
relators of negative sign, whereas a2^ is connected with 
positive further- distant correlation functions. Therefore, 
the difference in the sign of the correlators may be the 
reason for the relevance of the difference between ai^ and 
a2fi- This is in contrast to the FM case, where all corre- 
lators have a positive sign, and the equality a2^ = ai^ 
is a good assumption. Accordingly, we put a2^ = a2 
(cf. Ref. '22'), and, as in the FM case, we take A^ = A. 
To determine the five parameters aiy, aij^, 012, A and 
C at T = 0, we have the sum rule (fTT|) . the isotropy 
condition (fT^ . and the LRO condition ujq = 0. As 
the two additional conditions for fixing the free param- 
eters, we assume A(0) to be equal to the FM value, i.e., 
A(0) = 2 — and adjust the ground-state energy u(0) 
to the expression given by the linear spin-wave theory 
(LSWT) u{0) = ULswrjO) = -5(5 -H)( 2J,| + J^) + 
f Eq VC^JW + - (2^ll7q + J± cosg,)2. At finite 
temperatures, besides Eqs. ([TT|) and ([12]), and ujq — 
(for T ^ Tjv), we take Eq. ^ with a||(r) replaced by 
ai||(r) and the analogous condition (cf. Refs . [iSl and 1221) 



and 



ra[T) = — = ra[0). 



ai||(T)-l 



(19) 
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III. RESULTS 

As described in Sec. HH the quantities of the RGM de- 
termining the thermodynamic properties have to be nu- 
merically calculated as solutions of a coupled system of 
nonlinear algebraic self-consistency equations. For exam- 
ple, considering the antiferromagnet at T ^ T/v, we have 
11 equations for Cioo, Cqoi, Cno, C200, Cioi, C002 [ap- 
pearing in Eqs. I©-© and calculated by Eq. P^ ]. 
cti±, ct2, A and C. To solve this system of equations, we 
use Broyden's method, '^^ which yields the solutions with 
a relative error of about 10^^ on the average. The mo- 
mentum integrals occurring in the self-consistency equa- 
tions are done by Gaussian integration. 



A. Two-dimensional S ^ 1 magnets 

To test the quality of the approximations made in the 
RGM, in particular the assumptions about the vertex pa- 
rameters introduced in the decouplings ([3]) and (H]), we 
consider some correlation functions and thermodynamic 
properties of 2D S ^ 1 magnets in comparison with ED 
and QMC data. To provide a better comparison of the 
RGM with ED results, we apply the RGM also to finite 
systems with periodic boundary conditions proceeding 
as in Ref. 113. In Table U our RGM and ED results for 
several correlation functions of the 2D antiferromagnet 
at T = 0, also obtained by the RGM for a TV = 4 x 4 
square lattice, are presented. Determining the parame- 
ters (see Sec. |ll]) for the finite system with N — 16, as 
an input we take the ground-state energy in the LSWT 
that is also evaluated for N = 16. Let us consider the 
NN correlator Cio(O) determining the ground-state en- 
ergy u(0) = 3J||Cio. The LSWT and ED results are in 
a good agreement (for S = 2 they differ by only 0.1%). 
This provides some justification for using the LSWT data 
for u(0) as an input also in the 3D AF case. Note that 
the LSWT input is of advantage as compared with the 
choice made in Ref. 22, where u(0) is composed approx- 
imately from ID and 2D energy contributions which is 
justified for J±/J\\ ^ 1 only. The further-distant corre- 
lators listed in Table [J and calculated by the RGM for 

= 16 agree remarkably well (with an average deviation 




FIG. 1: 2D S — 1 ferromagnet: NN correlation function Cio, 
uniform static susceptibility x (upper inset), and specific heat 
Cv (lower inset) as functions of f = r/[| Jy \S{S + 1)], where 
the results of the RGM in the thermodynamic limit (solid 
lines) and for N = 8 (dashed lines) and the ED data (o, 
= 8) are shown. 



of 0.2%) with the ED results. 

Considering the 2D S = 1 ferromagnet, in Fig. [1] the 
temperature dependence of Cio, X; ai^d is plotted. 
For the finite lattice with A^ = 8, a very good agreement 
of the RGM and ED data is found. The comparison with 
the RGM results for N 00 demonstrates the finite-size 
effects. 

Next, we consider the 2D antiferromagnet at finite 
temperatures. Since the case S = 1/2 was intensively 
studied by the RGM in previous work )^^i^° we compare 
our results for S" = 1 with available QMC data.— As can 
be seen in Fig. [21 we obtain a surprisingly good agree- 
ment of the RGM with the QMC results (note that the 
QMC data for the correlation length agree with the SE 
results of Ref. Issh . This agreement is much better than 
for the S" = 1/2 antiferromagnetiii^ Correspondingly, for 
S = 1 we can give a rather reliable value for the zero- 
temperature susceptibility, x(0) = 0.07197. 

As outlined in Sec. [Ill i^i our approach more vertex 
parameters are introduced as independent equations for 



TABLE I: Correlation functions Cr of the 2D antiferromagnet at T = 0, as obtained by the RGM in the thermodynamic limit 
and for a finite system with A'' = 16, denoted by RGM(16), in comparison with the ED data for A*" = 16. 







S=l 






S=3/2 






S=2 




R 


RGM 


RGM(16) 


ED 


RGM 


RGM(16) 


ED 


RGM 


RGM(16) 


ED 


(1,0) 


-0.7720 


-0.7947 


-0.7980 


-1.6579 


-1.6920 


-1.6954 


-2.8773 


-2.9227 


-2.9261 


(1,1) 


0.5985 


0.6156 


0.6169 


1.3977 


1.4230 


1.4242 


2.5303 


2.5638 


2.5650 


(2,1) 


-0.5406 


-0.6032 


-0.6029 


-1.3109 


-1.4040 


-1.4035 


-2.4146 


-2.5383 


-2.5376 


(2,2) 


0.5077 


0.5649 


0.5689 


1.2616 


1.3462 


1.3503 


2.3488 


2.4611 


2.4651 
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FIG. 2: 2D S — 1 antiferromagnet: Correlation length ^, 
uniform static susceptibility x (upper inset), and staggered 
structure factor Cq = |Cq (lower inset) as functions of 
f = r/[J||S(S'+ 1)], where the RGM resuhs (solid lines) 
are compared with the QMC data of Ref. [s^ (•)• For com- 
parison, the results of a simplified version of the RGM with 
A(T) — A(0) (see text) are depicted (dot-dashed lines). 



them can be provided by the RGM. Therefore, we have to 
formulate appropriate additional conditions for their de- 
termination. Let us discuss, in comparison to the choice 
fixed in SecHIl two alternate choices of the parameters a2 
and A for the 2D S = 1 antiferromagnet (in two dimen- 
sions we omit the index /i =||, e.g., ai^2|| = «i,2), which 
are analogous to the choices made previously for the 
5=1/2 antiferromagnet and the S ^ 1 ferromagnetJ^ 
(i) If we choose a2 = ai, the parameter A(0) can be 
calculated (note that a2 and A only appear in the com- 
bination given by A) and used in Eq. (|18p . Then, we 
find the finite-temperature results to be not in such a 
good agreement with the QMC data as the results ob- 
tained by the parameter choice with a2 ^ ai- This corre- 
sponds to the findings for the S = 1/2 antiferromagnei^i^ 
and may be understood as explained in Sec. |TT1 There- 
fore, we discard the choice a2 — cti- (ii) If we adopt 
ct2 oil, but neglect the temperature dependence of 
A, i.e., A(T) = A(0) = 2 — -| (as was assumed for the 
FM case in Ref. '197, the results appreciably deviate from 
the QMC data, as is demonstrated in Fig. [2] (dot-dashed 
lines). This gives strong arguments for taking into ac- 
count the decrease of A(T) with increasing temperature 
[e.g., for 5* = 1, we have A(0) = 1 and Aoo = 0.625] and 
for our choice of the parameters for the antiferromagnet 
outlined on Sec.|TTl Note that for the 5* = 1 ferromagnet, 
where 0:2 = oli, the results shown in Fig. [T] only slightly 
improve those obtained by the assumption A(r) = A(0). 



B. Transition temperatures 

An important problem in the study of layered ferro- 
magnets and antiferromagnets is the calculation of the 
transition temperature Tm (M = C, TV) as a function 
of the interlayer coupling Jj_ and of the spin quantum 
number S. From the experimental side, the knowledge 
of the dependence Tm{R,S) with R = J±/J\\ is use- 
ful to estimate the interlayer exchange coupling from 
measurements of Tm- To test the quality of analyti- 
cal approaches, the precise results of numerical meth- 
ods, such as the QMC- and SE data^^ should be used 
as benchmarks. Considering the 3D isotropic model 
{R = 1), we have the inequahty^ Tn > Tc- More- 
over, Tm = Tm /[\J\\\S{S + 1)] is found to increase with 
increasing values of Considering, for example, the 
RPA, those results are not reproduced, instead we have 
fRPA ^ j.RPA^ .^i^gj.g J.RPA ig independent of S'-ii For 

layered magnets with R < 1, QMC and SE data in the 
FM case are still missing, so that there are no precise 
statements about the relation between Tjy and Tc as 
function of the interlayer coupling. With respect to the 
agreement with the QMC and SE data, our approach 
represents an important improvement as compared, e.g., 
to the RPA, which is outlined in the following. 

For the 3D ferro- and antiferromagnets, the solution of 
the RGM self-consistency equations yields the magneti- 
zation m(T) with m(TM) = at the second-order phase 
transition temperature Tm, where limj^^o Tm = is 
in agreement with the Mermin- Wagner theorem.'^'* In 
Fig. [3] and Table |TT] our results for Tm as functions of 



1.5 




FIG. 3: Neel temperature Tn = rjv/[J||S'(5'-|-l)] as a function 
of the interlayer coupling R = J±/J\\ ■ The results of the RGM 
(solid lines) and of the empirical formula (|20|l (dashed lines) 
are compared with the QMC data (•, Ref.|6i), the RPA (dotted 
lines, Ref. [Tib , and, for S = 1/2 (inset), with the mean-field 
theories of Refs.|ll (D),!!! (dot-dashed line), and[ll (dot-dot 
dashed line). 
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TABLE II: Transition temperatures Tm — Tm /[\J^^\S{S + 1)] of the ferromagnet (Tc) and antiferromagnet (Tjv) calculated by 
the RGM for different spins S and interlayer couplings J±/,J\\. 







Ferromagnet 




Antiferromaf 


^net 






5 = 1/2 


5 = 1 


S = 3/2 


C I/O 


C 1 
o = 1 


5 = 3/2 


S — oo 


0.0001 


0.2457 


0.3243 


0.3542 


0.1589 


0.3393 


0.3681 


3305 


0005 


2Q28 


0.3758 


0.4041 


0.2150 


0.4014 


0.4170 


0.3785 


0.001 


0.3184 


0.4027 


0.4298 


0.2498 


0.4331 


0.4421 


0.4039 


0.005 


0.3961 


0.4803 


0.5035 


0.3694 


0.5195 


0.5133 


0.4784 


0.01 


0.4403 


0.5226 


0.5436 


0.4430 


0.5640 


0.5521 


0.5200 


0.02 


0.4935 


0.5725 


0.5911 


0.5311 


0.6150 


0.5986 


0.5698 


0.05 


0.5826 


0.6552 


0.6706 


0.6681 


0.6979 


0.6776 


0.6538 


0.1 


0.6699 


0.7368 


0.7503 


0.7870 


0.7800 


0.7583 


0.7378 


0.5 


0.9953 


1.0571 


1.0694 


1.1655 


1.1121 


1.0884 


1.0667 


1.0 


1.2346 


1.3063 


1.3208 


1.4382 


1.3762 


1.3478 


1.3189 



R and S are presented, where in Fig. [3] the Neel tem- 
perature Tjv is compared with the QMC data of Ref. 6 
and other approaches. For the S ~ 1 antiferromagnet 
we get a very good agreement with the QMC results, 
as was also found for the 2D model (see Fig. [2]). Re- 
markably, the RPA results for both the 5' = 1 (Ref. [ll|) 
and S — 1/2 model o^°i^^ are in a rather good agreement 
with the QMC data. Considering the case 5* = 1/2 (in- 
set of Fig. [3]) and R < 0.04, we ascribe the reduction 
of Tn found by the RGM as compared to the RPA and 
the mean-field approaches of Refs. [H and [l^ to an im- 
proved description of strong AF quantum fluctuations at 
low temperatures counteracting the formation of LRO. 
For further comparison, the Neel temperature given very 
recently^^ by the interlayer mean-field approach within 
the Schwinger-boson mean-field theory is depicted for 
S — 1/2. The marked difference to the other curves 
(also found for S = 1) might be due to the asymmetry 
between intralayer and interlayer correlations introduced 
in this approach. 

Next we consider the transition temperatures Tj\/ for 
arbitrary values of S. The RGM yields fc{S) ^ TNiS), 
as can be seen in Table [TTl which is in accord with the 
QMC and SE data, but in contrast to the RPA result (see 
above). In passing to the classical limit S* — > oo we find 
limg^oo Tm — T^/"^ for all values of R. This may be un- 
derstood as follows. The RGM is a second-order theory 
that goes one step beyond the RPA and, therefore, pro- 
vides a better description of quantum fluctuations. Their 
vanishing for S oo may be reflected in the equality of 
the transition temperatures. 

We compare our results for the 3D isotropic model 
{R = 1) with the SE^ and QMC data^ for different spins. 
For the ferromagnet, the Curie temperatures Tc devi- 
ate from the SE values,^ Tc = 1.119 (1.2994, 1.37) for 
S" = i, (1, I), by 10% (0.5%, 4%). For the antifer- 
romagnet, the deviations of the Neel temperatures T/v 



from the SE values^ [agreeing with the QMC values for 
= 1/2 and = 1 (Ref. 6)], fjv = 1.259 (1.3676, 1.404) 
for 5" = i, (1, I), amount to 14%, (0.6%, 4%). From the 
experimental point of view, for the fit of exchange cou- 
pling parameters, deviations in the magnitude of tran- 
sition temperatures of up to about 10% are considered 
as a reasonable accuracy. In both the FM and AF cases 
the RGM yields the best values of Tm for 5 = 1. For 
any spin, we get the correct relation Tjy > Tc, where 
the ratio Q = T^/Tc = 1.17 (1.05, 1.02) for 5=1, 
(1, |) agrees well with the SE values Q = 1.13 (1.05, 
1.03). That means, concerning the difference between 
Tjv and Tc, the RGM yields good results for all values 
of S. Considering the dependences Tm{S), the increase 
of Tc with increasing S is in qualitative agreement with 
the SE data. For the antiferromagnet, TV decreases with 
increasing S being opposite to the behavior of the SEi^ 
and QMC dataj^ This is connected with the inequality 
Tn^S = I) > lims^ooTW = 1.3189, whereas the QMC 
data&^ yield T'Ar(5 = i) < fjv(5 = oo) = 1.443 [note 
that in the classical Hciscnberg modeP^ the spins are 
taken of unit length, and the exchange interaction J'^' is 
related to J = = J,, by = JS{S + 1)]. 

Let us consider the anisotropic magnets {R < 1). For 
5 = 1/2 and i? < 0.01 we find T^ <Tc, and for R ^ 0.01 
we have T/v > Tc- In the cases S — 1 and S — 3/2 we get 
T/v > Tc for all values of R. The peculiarity in the rela- 
tion between T/v and Tc for 3=1/2 may be explained 
by the presence of strong AF quantum fluctuations at 
low temperatures which may suppress the AF LRO. 

For the discussion of experimental data it is convenient 
to use an analytical expression for Tm{R, S). Our RGM 
results for the dependence of Tm on R may be well fitted 
by the empirical formula proposed in Ref. |^, 



Tm = 



A 



B-\n{J^/J\\y 



(20) 



7 



TABLE III: Coefficients of the empirical law [Eq. pOfl ] for the transition temperatures of the ferro- and antiferromagnet. 







Ferromagnet 




Antiferromagnet 






5=1/2 


S = 1 


5 = 3/2 


5 = 1/2 5 = 1 5 = 3/2 


5 = oo 


A 


3.15 


4.00 


4.27 


1.95 4.36 4.34 


3.96 


B 


2.50 


3.08 


3.27 


0.01 3.21 3.27 


3.01 



where the values of A and B are listed in Table IIIII The 
concrete values of the coefficients slightly depend on the 
choice of data points used for the fit. Since Tm reveals 
the strongest increase with R for i? <^C 1 , in this region we 
take points lying more dense than for moderate interlayer 
couplings. The values given in Table Hill are obtained by 
choosing points within the interval R = 10"'' to 10"^ and 
R = 10-2 to 1 in steps of Ai? = lO'^ and Ai? = 10"^ 
respectively. Then, a good fit in the whole R region can 
be achieved in all cases, except for the 5=1/2 antifer- 
romagnet, where a reasonable fit by Eq. (PO)) is obtained 
for i? ^ 0.1 (see Fig. [3]). 



C. Specific heat 

The temperature dependence of the specific heat Cv 
is characterized by a cusplike singularity at the transi- 
tion temperature Tm determined by and, for suffi- 
ciently low interlayer couplings, by a broad maximum 
above Tm that is mainly determined by Jn . For the 3D 
isotropic magnets, Cv is plotted in Fig. HI Considering 
the S = \/2 ferromagnet (see inset), above Tc we ob- 



tain an excellent agreement with the SE data of Ref. 
For the 5 = 1/2 antiferromagnet, the agreement of the 
RGM with the QMC results^ is very good at tempera- 
tures sufficiently below and above T/v, whereas near T^r 
the height of the cusp is underestimated. Considering 
the S dependence of Cy in the LRO phase, with increas- 
ing S the slope of the Cv curves near Tm decreases, and 
the cusp develops to a kink (see Fig. [5]) . The analogous 
tendency is found in the FM case. This behavior may 
be considered as a deficiency of the RGM, because in the 
classical Heisenberg model {S — > oo) the QMC data of 
Ref. sd yield evidence for a cusplike structure of Cv at 
Tm- 

Next we consider the specific heat of quasi-2D magnets. 
In the ferromagnet a broad maximum, in addition to the 
phase-transition singularity, appears at i? < 0.035 (i? < 
0.015) for S = 1/2 {S = 1), as can be seen in Fig.lH The 
analogous behavior is found for the antiferromagnet, as 
shown in Fig.[Sl Here, the broad maximum occurs at i? < 
2-3 (i? <; 0.015) for S = 1/2 (S=l), which agrees with 
the 5=1/2 QMC data of Ref. i. As for the isotropic 
5=1/2 antiferromagnet (cf. Fig. 2]), the RGM agrees 
well with the QMC results at low and high temperatures. 




0.35 



d 0.25 h 




0.15 



T/T, 



FIG. 4: Specific heat Cv of the isotropic antiferromagnet 
(Jx = J||) with 5 = 1/2, 1, 2, and 5, from bottom to top, 
in comparison to the QMC data for 5 = 1/2 (•, Ref. |7|). 
The inset displays Cv for the isotropic 5=1/2 ferromagnet, 
compared to the SE results of Ref. i2l (♦)• 



FIG. 5: Specific heat Cv of the ferromagnet with 5=1/2 and 
5 = 1 (inset) for Jx/Jy = 0.01 (solid lines), with 5 = 1/2 for 
Jx/J|l = 0.025 (dot-dot dashed line) and 0.035 (dotted line), 
and with 5 = 1 for Jx/Jy = 0.015 (dashed line) and 0.02 
(dot-dashed line). 
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1.00 - S=l/2 



0.75 - 



0.50 - 




0.25 - ; 



o.oo„ 



FIG. 6; Specific heat of tlie 5=1/2 antiferromagnet in 
comparison to tlie QMC data of Ref. 5: (filled symbols) for 
Jx/J|| = 2"^ (solid line, •), 2"^ (dashed line, ►), and 2~^ 
(dot-dashed line, The inset shows Cv of the S = 1 anti- 
ferromagnet for J±/ J\\ = 0.01 (solid line), 0.015 (dashed line), 
and 0.02 (dot-dashed line). 



Again, the height of the cusp is underestimated, where 
the relative deviation of Cv{Tn) from the QMC values 
increases with decreasing R. 

Recently, specific heat data for the quasi-2D 5=1/2 
antiferromagnet Zn2VO(P04)2 were presentedi^ Taking 
Tn = 3.75K and J,, = 7.41K from Ref. M, by Eq. ^ 
and Table Hill we get i? = 5.8 x lO"^. Calculating the 
specific heat we obtain a broad maximum at Tm = 5.9K 
with Cv{Tm) — 0.45 which corresponds to the measured 
broad hump at Th = 4.5K with the height Cv{Th) = 0.45 
agreeing with the theoretical value of Cv{Tm)- At T/v, 
the experiment shows a pronounced cusp with Cv{Tn) — 
0.6. As discussed above (see Fig. [S]), this feature cannot 
be reproduced by the ROM, instead we get a small spike 
at Tn with Cv{Tn) ^ 0.3. 



D. Correlation length 

The intralayer and interlayer correlation lengths 
{fi =11, ±) for R ^ diverge as T approaches Tm from 
above. In the vicinity of Tm, S,^^ and behave as T — 
Tm (corresponding to the critical index v = 1) also found 
by previous mean-field approaches J^ii^ This can be seen 
in Fig. [3 that shows versus f = T/[\J\i\S{S + 1)] of 
the 5=1/2 and 5=1 ferromagnet. The curves for the 
antiferromagnet look similar. At fixed R < 1 and 5 we 
have < which corresponds to the weaker interlayer 
as compared to the intralayer correlations. Considering 
the 5 dependence of for the ferromagnet, we have 
fb(5 = 1/2) < fciS = 1) (see Fig. [7] and Sec. [IIlB|) 
which implies, at fixed T > Tc and R, the inequality 



1.5 



1 - 



0.5 





[ 1 1 y 








' / 




I /7 






























II.// 1 . 




FIG. 7: Inverse correlation lengths within (^y ^) and between 

the xy planes (Cx"^, see inset) versus T = r/[| J|| |S'(S' + 1)] of 
the ferromagnet with 5 = 1/2 (solid lines) and 5 = 1 (dashed 
lines) for J±/J|| = 0.001, 0.1, and 1, from left to right. 



^^(5 = 1/2) < ^^(5 = 1). Note that recently, an analo- 
gous 5 dependence for the lonitudinal correlation length 
of the 2D ferromagnet in a small magnetic field was 
found, also by QMC,.2^ i.e., ^''''[S = 1/2) < ^^(5" = 1) 
at fixed f. 

Let us compare our results for the intralayer corre- 
lation length ^[| with the neutron-scattering data on 
the 5=1 quasi-2D antiferromagnet La2Ni04i^ Tak- 
ing Tn = 327.5K and Jii = 28.7meV from Ref. i, by 




FIG. 8: Antiferromagnetic intralayer correlation length in 
La2Ni04 obtained by the neutron-scattering experiments of 
Ref. (o) compared to the QMC data for Jx = (•, 
Ref. 1321 ) and the RGM results for Jx = (solid line) and 
Jx/J|| = 3.5 X 10"^ (dashed line). 
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Eq. ^ and Table |TTT] we obtain i? = 3.5 x IQ-^. In 
Fig. [5] the experimental data are plotted in comparison 
to the QMC data for i? = (Ref. 'M) and the RGM 
results for R = and R = 3.5 x 10"'^, where a sat- 
isfactory overall agreement with experiments is found. 
At fixed temperature, the correlation length for R > 
is larger than for R = 0, because ^|| diverges at T/v. 
To explain the neutron-scattering experiments, in Ref. |^ 
a small Ising anisotropy in the strictly 2D model was 
considered which leads to in a finite transition tempera- 
ture somewhat below Tpf. Such an easy- axis anisotropy 
was also discussed in Ref. ITtI to explain the experiments. 
However, as was shown in Ref. [s^, the experimental data 
with ^"''P < ^QMC (|ggg pjg ^ g^j.g incompatible with the 
QMC results obtained for the 2D model with a small Ising 
anisotropy, since it even enhances the correlation length 
at low temperature. In our approach, the finite value of 
Tjv is ascribed entirely to the interlayer coupling which 
gives e > ^QMC^ To improve the agreement with exper- 
iments, let us point out, that in our calculations a sim- 
ple cubic lattice was taken, whereas in the orthorhombic 
structure of La2Ni04 the interlayer coupling is frustrated. 
As was shown in Ref. [29|, in the Ji — J2 model, frustration 
may appreciably reduce the correlation length. The in- 
fluence of frustration on the transition temperature and 
correlation length of quasi-2D Heisenberg magnets will 



be left for further study. 



IV. SUMMARY 

In this paper the thermodynamics of layered Heisen- 
berg magnets with arbitrary spin S is systematically 
investigated by a spin-rotation-invariant Green-function 
method and by exact diagonalizations on finite 2D lat- 
tices. The main focus is put on the calculation of the 
Curie temperature Tc and the Neel temperature Tjv 
in dependence on the interlayer coupling J± and the 
spin quantum number. From the numerical data we ob- 
tain simple empirical formulas for Tc.n{J±)- A good 
agreement of our results, in particular on the relation 
between Tc and T^r, with available quantum Monte 
Carlo and series-expansion data is found. The compar- 
ison to experiments on the quasi-2D antiferromagnets 
Zn2VO(P04)2 and La2Ni04 yields a reasonable agree- 
ment. From our results we conclude that the application 
of the second-order Green-function approach to extended 
layered Heisenberg models (frustration, anisotropy in 
spin space) may be promising to describe the unconven- 
tional magnetic properties of real low-dimensional quan- 
tum spin systems. 
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